DISTRIBUTION  STATEMENT  A.  Approved  for  public  release;  distribution  is  unlimited. 


Title:  4DVar  NCOM  in  RELO  and  COAMPS5 

PI:  Dr.  Hans  Ngodock 
1009  Balch  Blvd 
Stennis  Space  Center,  MS  39529 

Tel:  (228)  688  5455  Fax:  (228)  688  4759  Email:  Hans.Ngodock@nrlssc.navy.mil 

Co-PI:  Dr.  Matthew  Carrier 
1009  Balch  Blvd 
Stennis  Space  Center,  MS  39529 

Tel:  (228)  813  4086  Fax:  (228)  688  4759  Email:  Matthew.Carrier@nrlssc.navy.mil 

Co-PI:  Dr.  Clark  Rowley 
1009  Balch  Blvd 
Stennis  Space  Center,  MS  39529 

Tel:  (228)  688  5809  Fax:  (228)  688  4759  Email:  Clark.Rowley@nrlssc.navy.mil 

Co-PI:  Dr.  Scott  Smith 
1009  Balch  Blvd 
Stennis  Space  Center,  MS  39529 

Tel:  (228)  688  4630  Fax:  (228)  688  4759  Email:  Scott.Smith@nrlssc.navy.mil 

Co-PI:  Dr.  Timothy  Campbell 
1009  Balch  Blvd 
Stennis  Space  Center,  MS  39529 

Tel:  (228)  688  5285  Fax:  (228)  688  4759  Email:  Scott.Smith@nrlssc.navy.mil 
Award  number:  N0001413WX20550 


LONG  TERM  GOALS 

To  improve  the  accuracy  of  the  ocean  data  assimilation  and  forecast  by  implementing  the  four¬ 
dimensional  variational  (4DVAR)  technique  with  the  Navy  coastal  ocean  model  (NCOM) 

OBJECTIVES 

This  project  will  implement  a  4DVar  capability  into  the  relocatable  /  regional  NCOM  system  and  the 
COAMPS5  system  as  well  as  merge  the  NCODA  and  NAVDAS-AR  frameworks.  The  tangent  linear 
and  adjoint  for  NCOM  will  be  for  the  latest  version  using  the  generalized  vertical  coordinate  system. 
All  code  will  be  added  to  the  appropriate  subversion  repositories.  The  code  will  be  tested  on  the 
operational  computer  systems  prior  to  a  validation  test  report.  The  results  will  be  validated  in  one  area 
to  be  determined  in  coordination  with  the  operational  POC.  The  area  extent  of  the  validation  area  will 
be  balanced  against  available  computational  resources  in  conjunction  with  the  operational  POC.  A 
validation  test  report  for  the  area  relative  to  available  in  situ  data  will  be  written.  The  evaluations  will 
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be  for  temperature  structure,  salinity  structure  and  acoustic  proxy  parameters  relative  to  in  situ  profile 
data.  In  addition,  evaluations  will  be  constructed  for  velocity  comparisons  to  in  situ  drifter  data.  All 
evaluations  will  be  constructed  with  the  new  4DVar  system  relative  to  the  existing  operational  3DVar 
system.  The  code  will  be  available  to  all  Navy  operational  centers  to  OPTEST. 

APPROACH 

The  approach  is  divided  into  the  following  tasks 

Integration  with  Regional  NCOM  and  COAMPS 

At  present,  the  4DVar  uses  the  output  files  generated  by  the  NCODA  OcnQC  data  prep.  The  4DVar 
must  be  integrated  with  the  merged  NCODA/NAVDAS-AR  code.  Initially  the  4DVar  will  be  called 
from  NCODA  so  that  tests  can  proceed  on  the  HPCMP  platfonns  and  other  work  can  progress.  As  the 
NCODA  and  NAVDAS-AR  are  merged,  the  4DVar  will  also  be  incorporated. 

Parallelization  and  scalability 

The  present  parallelization  of  the  4dvar  code  leverages  the  parallelization  of  the  NCOM  code  itself.  It 
uses  the  same  tiling  for  the  TL  and  adjoint  models  as  in  the  parent  NCOM  code.  The  conjugate 
gradient  step  is  perfonned  simultaneously  by  all  processors. 

Sigma-z  and  generalized  vertical  coordinates 

The  4dvar  code  was  developed  with  the  premise  that  the  generalized  vertical  coordinate  (gvc)  of 
NCOM  was  going  to  be  the  standard  operational  code.  The  transition  of  the  gvc  NCOM  to  NAVO  is 
under  way  presently.  Given  that  NCOM  can  run  in  gvc  mode  with  a  sigma-z  setup  and  vice  versa,  we 
need  to  ensure  that  this  property  is  maintained  in  the  4dvar  code,  i.e.  the  gvc-based  4dvar  can  run  with 
a  sigma-z  setup. 

Strong  or  weak  constraints 

In  an  effort  to  lower  the  computational  cost  of  the  4dvar,  a  strong  constraint  approach  can  be 
implemented.  A  pure  strong  constraint  form  will  need  a  steepest  descent  or  quasi-Newton 
minimization  not  yet  implemented,  whereas  in  the  non-pure  form  we  just  set  the  model  error  term  to 
zero  in  the  system  we  have  now.  We  will  investigate  the  difference  in  both  approaches  and  quantify  the 
extra  computational  cost  incurred  by  the  weak  constraint  approach. 

Cost  effectiveness  to  meet  time  constraints 

Can  the  4dvar  analysis  be  done  in  a  specific  amount  of  time,  say  2  hours?  This  will  depend  on  the 
region,  the  preconditioning,  and  mostly  the  time  window.  Balances  of  these  parameters  will  be 
constructed  to  ensure  that  the  system  fits  within  operational  constraints. 

Covariance  space  and  relative  errors  of  observations  and  initial  conditions 

The  horizontal,  vertical  and  temporal  scales  of  the  covariances  have  to  be  prescribed  carefully 
according  to  the  prevailing  dynamics  over  the  chosen  region.  The  covariance  model  in  place  is  based 
on  the  solution  proposed  by  Egbert  et  ah,  (1994)  and  Weaver  and  Courtier  (2001)  using  a  spatial 
transfonnation  of  a  Gaussian  function  to  represent  spatial  structure  dependent  on  any  parameter 
desired.  This  was  extended  to  use  the  diffusion  solution  (Bennett  et  al.,  1996,  1998,  Ngodock  2005)  to 
apply  the  covariance  efficiently  (Carrier  and  Ngodock,  2010,  Yaremchuk  et  al.,  2011).  At  present  the 
vertical  structure  is  uncorrelated,  and  the  vertical  correlation  structure  must  be  added.  Also,  as  we 
have  found  in  the  Monterey  Bay,  the  correlation  scales  have  a  direct  impact  on  the  assimilation  and 
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forecast  skills.  This  task  consists  of  selecting  the  correlation  scales  as  well  as  implementing  the  vertical 
correlation  model  and  evaluating  them  in  tenns  of  analysis  and  forecast  skills. 

Assimilation  cycle 

The  4dvar  system  is  expected  to  fit  the  data  within  the  time  range  of  TLM  stability.  The  latter  depends 
on  resolution,  dynamical  regime  and  geographic  location.  Given  the  computational  demands  of  the 
system,  a  short  assimilation  cycle  will  most  likely  fit  the  time  allocated  to  perform  the  analysis.  There 
has  not  been  a  concise  study  to  determine  what  the  length  of  the  assimilation  should  be  depending  on 
the  considerations  enumerated.  There  are  several  approaches  that  need  to  be  understood  to  properly 
implement  this.  First,  it  is  possible  to  use  several  1-day  cycles  in  hind  cast  versus  a  single  multiday 
hind  cast  cycle,  ft  is  also  possible  extend  the  analysis  period  into  the  future  so  that  observation  could 
be  projected  into  corrections  beyond  the  final  observation  time.  Evaluations  of  these  alternative 
implementations  will  be  made,  and  the  operational  implementation  will  be  based  on  the  results  of  the 
tests. 

Demonstration  that  the  4DVAR  improves  the  forecast  skill  of  3DVAR 

There  needs  to  be  a  robust  demonstration  that  the  4dvar  is  an  improvement  to  the  3dvar  and  that  it 
justifies  the  computational  cost.  This  demonstration  has  to  be  done  in  the  area  agreed  upon  with  the 
operational  POC.  The  evaluations  relative  to  present  capabilities  in  3DVar  will  be  part  of  the 
validation  test  report. 

WORK  COMPLETED 

During  FY13,  the  4DVar  has  been  successfully  integrated  into  the  relocatable  NCOM  (RELO).  This 
work  involved  modifying  the  RELO  scripts  to  allow  for  running  the  4DVar  analysis  suite  in  place  of 
the  3DVar.  In  addition  to  this,  minor  modifications  were  made  to  the  NCODA  OcnQC  data  prep 
programs  to  cycle  observations  for  a  time-window.  Work  to  integrate  4DVar  into  COAMPS-5  is 
ongoing  with  expected  completion  by  the  end  of  QTR1  FY14. 

The  4DVar  has  been  proven  to  work  for  any  vertical  grid  that  can  be  specified  for  NCOM,  including 
the  standard  sigma-z  and  the  generalized  vertical  coordinate.  Earlier  issues  with  TLM  stability 
regarding  deep  sigma  layers  were  resolved  in  QTR2  FY13  by  a  minor  modification  to  the  linearization 
of  the  vertical  mixing  parameters.  The  4DVar  is  now  stable  for  any  NCOM  configuration  and 
resolution. 

As  of  FY13,  the  4DVar  is  fully  parallelized.  Tests  of  the  4DVar  on  the  DSRC  machines  (KILRA1N 
and  LIAISE)  show  that  the  scalability  of  the  4DVar  matches  that  of  NCOM  (albeit  at  a  higher 
computational  cost  overall).  Work  is  continuing  into  FY14  to  determine  the  best  operational 
configuration  of  the  4DVar  that  is  a  balance  between  the  quality  and  computational  efficiency  of  the 
analysis.  This  work  includes  investigating  the  assimilation  cycling  (i.e.  length  of  the  data  window)  as 
well  as  the  errors  for  the  initial  condition  and/or  model  fields.  All  of  these  aspects  affect  the  quality 
and  efficiency  of  the  analysis,  and  therefore  must  be  investigated  in  this  light.  Expected  completion  of 
this  work  is  the  middle  of  QTR1  FY14. 

In  addition  to  the  work  done  on  the  cycling  and  relative  error  settings,  a  preconditioner  has  been  added 
to  the  minimization  software  used  in  the  4DVar.  This  preconditioner  (based  on  Courtier  1997)  has 
resulted  in  a  significant  cost  savings  by  reducing  the  number  of  minimization  steps  needed  to  achieve 
convergence  (figure  1). 
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The  assimilation  of  sea  surface  height  (SSH)  has  proven  problematic  in  shallow  water  regions.  In 
these  locations  (depths  less  than  200m),  the  assimilation  of  SSH  results  in  the  excitation  of  gravity 
waves,  which  degrade  the  resulting  analysis.  We  have  temporarily  addressed  this  issue  by  increasing 
the  spatial  correlation  for  the  SSH  field  in  the  assimilation,  thereby  smoothing  the  increment  added  to 
the  tangent  linear  model.  In  the  long-run,  however,  filtering  the  gravity  waves  in  the  solution  of  the 
free-surface  may  be  best;  this  will  be  incorporated  into  the  tangent  linear  and  adjoint  models  in  QTR1 
FY14. 

We  have  had  great  success  in  assimilating  velocity  observations  with  the  4DVar.  Velocity 
observations,  by  nature,  are  highly  variable  in  space  and  time.  This  aspect  makes  them  very  difficult  to 
assimilate  via  a  3DVar  method,  which  assumes  all  the  data  were  taken  at  one  time.  For  the  recent 
Grand  Lagrangian  Deployment  (GLAD)  experiment  that  took  place  in  the  Gulf  of  Mexico  in  the 
summer  of  2012,  a  large  number  of  surface  drifters  were  released.  Eulerian  velocity  values  have  been 
derived  from  the  position  data.  These  derived  velocity  observations  are  assimilated  into  the  4DVar 
(using  a  6  km  Gulf  of  Mexico  domain,  where  the  observations  are  assimilated  every  hour  during  a  4- 
day  data  window).  The  resulting  forecast  shows  excellent  skill  in  the  velocity  fields  out  to  96-hours 
(when  compared  to  a  forecast  that  did  not  assimilate  the  velocity  data).  Figure  2  shows  the  comparison 
between  selected  modeled  (purple)  and  observed  (green)  drifter  tracks  between  15-19  August,  2012  in 
the  northern  Gulf  of  Mexico.  The  top  (bottom)  panel  of  fig.  2  shows  the  comparison  using  an  NCOM 
forecast  that  excluded  (included)  GLAD  velocity  observations  in  the  analysis  step.  The  forecast  that 
used  the  analysis  that  included  the  velocity  observations  matches  the  actual  drifter  positions  very 
closely,  clearly  indicating  that  the  assimilation  of  velocity  observations  improves  the  prediction  of  the 
velocity  field. 

Early  comparisons  between  the  4DVar  and  3DVar  have  been  made  in  the  Okinawa  trough  domain 
selected  for  this  project.  Initial  indications  are  that  the  4DVar  produces  a  forecast  that  is  more  accurate 
than  the  3DVar,  in  terms  of  model-observation  differences  during  the  forecast  cycle  of  1-30  August, 
2007.  Figure  3  shows  the  absolute  difference  between  model  and  observed  temperature  for  no 
assimilation  (top  panel),  3DVar  assimilation  (middle  panel),  and  4DVar  assimilation  (bottom  panel); 
figure  4  shows  this  same  comparison  for  salinity.  All  the  model-observation  differences  for  figs.  3  and 
4  are  shown  for  profile  observations,  where  the  profiles  have  been  arranged  chronologically  (not 
spatially).  Fig.  3  shows  that  the  forecast  using  the  4DVar  analysis  has  lower  error  in  the  100-400  m 
depth  range  than  its  3DVar  counterpart.  A  similar  conclusion  can  be  made  for  salinity  in  fig.  4,  as  the 
error  here  is  lower  for  the  4DVar  than  the  3DVar.  It  is  important  to  note  that  these  results  are 
preliminary;  a  more  extensive  comparison  between  the  3DVar  and  4DVar  will  be  made  as  part  of  the 
validation  test  report  (VTR);  experiments  for  the  VTR  are  on-going  with  expected  completion  QTR  2 
FY14. 
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RESULTS 


Conjugate  gradient  convergence  OKT  1-6  Aug.  2007 


-  No  preconditioner  CG  iterations 

—  —  Preconditioner 


Conjugate  gradient  convergence  OKT  6-11  Aug.  2007 


Figure  1:  Conjugate  gradient  convergence  comparison  between  the  non-preconditioned  (solid  line) 
and  preconditioned  (dashed  line)  NCOM-4DVAR  for  5-day  assimilation  window  between  1-6 
August,  2007  (top  panel)  and  6-11  August,  2007  (bottom  panel). 
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Figure  2:  Comparison  of  modeled  (purple)  and  observed  (green)  drifter  tracks  in  the  northern  Gulf 
of  Mexico  between  15-19  August,  2012.  The  top  (bottom)  panel  shows  this  comparison  for  modeled 
drifters  using  NCOM forecast  velocity  that  excluded  (included)  velocity  measurements  in  the 
analysis  step.  Each  plot  is  overlaid  on  the  respective  SSH  field  for  each  forecast. 
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Data  minus  NCOM  (No  Assimilation)  Temperature  (deg  C) 


13  21  24  25 


Data  minus  NCOM  (3dvar  Assimilation)  Temperature  (deg  C) 


Data  minus  NCOM  (4dvar  Assimilation)  Temperature  (deg  C) 


13  21  24  25 

Days  of  August,  2007 


Figure  3:  Absolute  difference  between  NCOM  and  observed  temperature  profiles  in  the  Okinawa 
Trough  domain.  Observed  minus  NCOM  with  (1)  no  assimilation  (top),  (2)  3DVar  assimilation 
(middle),  and  (3)  4DVar  assimilation  (bottom)  is  shown.  All  profiles  are  plotted  chronologically  (not 

spatially)  from  1-30  August,  2007. 
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Data  minus  NCOM  (No  Assimilation)  Salinity  (psu) 


04  13  20  25 

Data  minus  NCOM  (First- Guess)  Salinity  (psu) 


04  13  20  25 

Days  of  August,  2007 


Figure  4:  Same  as  in  fig.  3,  but  for  salinity. 

IMPACT/APPLICATIONS 
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